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Different interaction models treating hadronic collisions below i?iab<100 GeV in air shower simulations are com- 
pared within the frame of the CORSIKA program. Their behaviour is studied in isolated collisions of protons or 
pions with light-mass targets in comparison with experimental data. It is discussed how model properties influence 
air-shower parameters like the longitudinal profile or the lateral distribution of particles arriving at ground. 



1. INTRODUCTION 

In extensive air showers (EAS) induced by cosmic 
hadronic particles low-energy secondary hadrons 
collide with the atmospheric nuclei thus forming 
the final branches of the hadronic shower back- 
bone. Recently Ref. ^U] reported that in sim- 
ulations the lateral particle distribution at core 
distances >1 km shows a strong dependence on 
the used low-energy hadronic interaction model. 
Therefore several codes describing low-energy in- 
teractions have been examined with the EAS sim- 
ulation program CORSIKA m the present 
study. Preliminary results are given in The 
used interaction codes are presented in Sect.EI In 
Sect.|3]their predictions are compared with exper- 
imental results of single hadronic collisions with 
light targets close in mass to the atmospheric con- 
stituents. How the behaviour in single collisions 
transforms into measurable properties of EAS is 
discussed in Sect. 0] The consequences to simu- 
lations are exemplified in Sect. 

2. MODELS 

Within this study the 4 codes FLUKA [HI, 
GHEISHA [H], Hillas Splitting Algorithm H3, 
and UrQMD |H] are investigated in view of 
their application to low-energy (i?i a b<100 GeV) 
hadronic interactions within the EAS simulation 
program CORSIKA. 

In some comparisons the high-energy interac- 
tion models DPMJET 11.55 EU, NeXus 3 [§], 
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QGSJET 01 HU, or SIBYLL 2.1 H2] are included 
overstretching their recommended energy range 
to show up possible problems in the transition 
region around 80 GeV. 

The hadronic event generator FLUKA 2002 [HI 
is used with CORSIKA only for the description 
of the inelastic interactions with laboratory ener- 
gies below several 100 GeV. Within FLUKA these 
collisions are handled by different hadronic inter- 
action models above, around, and below the nu- 
clear resonance energy range. The capability of 
FLUKA for cosmic ray calculations has recently 
been demonstrated 0. 

The GHEISHA ^H) program successfully used in 
the detector Monte Carlo code GEANT3 HHj is 
called GHEISHA 600 within this study to distin- 
guish it from version GHEISHA 2002 with modi- 
fied kinematics using correction patches |HJ which 
improve energy and momentum conservation. 
The Hillas Splitting Algorithm (HSA) EJJ is em- 
ployed in the EAS simulation code AIRES |26I27| . 
Its parametrizations are valid only for collisions 
with air which do not allow a direct comparison 
with experiments performed with other targets. 
The HSA has been linked with CORSIKA for test 
purposes only to add the results of these tests to 
the present comparison. It is not planned to make 
the HSA available within CORSIKA. 
The Ultra-relativistic Quantum Molecular 
Dynamic (UrQMD 1.3) model |H] describes mi- 
croscopically the projectile transport through the 
air target in tiny steps (sa 0.2 fm) and follows 
collisions and/or scatterings on the hadron level. 
Nuclear resonance effects are considered in detail. 
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Figure 1. 7r + (upper panel) and 7r~ (lower panel) 
multiplicities in pp-collisions as function of pro- 
jectile energy. Experimental points are taken 
from JSESEHl. 



3. SINGLE COLLISIONS 

All simulations presented in this section have 
been performed with the so-called 'interaction- 
test' option of CORSIKA to compare all interac- 
tion programs under equal conditions. 

3.1. pp-collisions 

The 7r + and 7r~ -multiplicities in pp-collisions are 
shown in Fig.^ Most of the displayed low-energy 
models follow the experimental points sufficiently 
well except GHEISHA which produces too few pi- 
ons. QGSJET overestimates the pion multiplicity 
within the 20 - 500 GeV region at the lower end 
of its design energy range. 

3.2. p- 9 Be collisions 

The above mentioned deficit of pions produced by 
GHEISHA is visible in the longitudinal momen- 
tum fractions of p- 9 Be collisions (Fig. |2J) at xi a b- 
values around 0.15. As in pp-collisions QGSJET 
shows a pion overestimation visible in Fig. at 
a^iab < 0.1. Despite the completely different ap- 
proaches of FLUKA and UrQMD both models 
reproduce the experimental data equally well. 
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Figure 2. Distributions of 7r + (upper panel) and 
7r~ (lower panel) longitudinal momentum frac- 
tions iriab = Ptot/Pbeam in p- 9 Be collisions at 24 
GeV/c. Experimental data were derived by Engel 
et al. [12 from |1I2U1 |. 



For kaons all codes behave similar as for the 
pions. The dominance of K + -mesons over the 
K~ ones caused by the associated production of 
K + +A in proton initiated interactions is repro- 
duced by all studied models similarly (shown in 
Ref. 19 1). FLUKA produces a small deficit of K ± 
at cci a b < 0.3. 

The agreement of the transverse momenta with 
experimental data is examined in Fig. where the 
double-differential cross sections of 7r + -mesons 
are given as function of transverse mass and ra- 
pidity. The rapidity bins of width Ay = 0.2 ex- 
tend from 0.6 (curve a at bottom) to 2.8 (curve 
k on top). Both FLUKA and UrQMD follow the 
exponential trend (straight lines in logarithmic 
scale) of the experimental points J] with the cor- 
rect slope in the whole rapidity range. GHEISHA 
exhibits significant deviations for small transverse 
masses mj_— mo < 0.2 and for rapidities 1.0 < y < 
2.4 (curves c to i). Presumably the kinematics are 
still far from being perfect despite the correction 
patches |S] . As the high-energy models QGSJET 
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Figure 3. Invariant cross sections as function of 
mj_ for 7T + in p- 9 Be collisions at 14.6 GeV/c. The 
curves a - k give the distributions for consecutive 
rapidity intervals Ay = 0.2 increasing in the range 
0.6<y<2.8. Each curve is multiplied by a factor 
10' with I increasing from to 10 to separate the 
curves. Experimental points are taken from pQ. 



and neXus use a p^-parametrization covering a 
wide energy range up to highest energies it is not 
surprising that below the lowest end of their de- 
signed energy range they produce transverse mo- 
menta with slopes generally too steep. All these 
findings hold similarly for 7r~-mesons emerging 
from p- 9 Be collisions. 

3.3. p- 12 C collisions 

For the carbon target the diagnoses resemble the 
conclusions from Sect. 13.21 Longitudinal momen- 
tum fractions of tt~ -mesons are plotted for pro- 
jectile energies of 4.2 and 10 GeV in Fig.0] Again 
the agreement of FLUKA and UrQMD with the 
experimental data is fine. Simultaneously for low- 
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Figure 4. Distribution of 7r~ longitudinal mo- 
mentum fractions a;iab = Ptot/Pbcam in p- 12 C col- 
lisions at 4.2 (upper panel) and 10 GeV/c (lower 
panel). Experimental data are taken from 4 . 



energy pions (cciab < 0.3) QGSJET produces an 
overshoot and GHEISHA shows a deficit at both 
examined energies. 

3.4. Collisions with 14 N and air targets 

In contrast to the preceding sections no exper- 
imental data exist for 14 N or air targets. We 
take the behaviour of FLUKA and UrQMD as 
reference as both models have shown a suffi- 
ciently good agreement with experimental values 
for lighter targets. 

As charged pions are by far the most frequent 
projectiles within an EAS and 14 N is the most 
abundant target nucleus within air, we investi- 
gate in Fig. the momentum fractions which the 
7r ± -mesons carry away from 7r + - 14 N collisions. 
The behaviour of the displayed models resem- 
bles that which has already been found for the 
lighter targets. FLUKA and UrQMD are close to 
each other. All other high energy models, at 20 
GeV collision energy neXus, at 100 GeV addi- 
tionally DPMJET, QGSJET, and SIBYLL result 
in similar xi a b-distributions (Again the overesti- 
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Figure 5. Distribution of ^ longitudinal mo- 
mentum fractions xi ab = Ptot/Pbcam in 7r+- 14 N 
collisions at 20 (upper panel) and 100 GeV (lower 
panel) . 



mation of QGSJET at low xiab-values deviates 
from the average behaviour). The peak occur- 
ring at xiab > 0.9 results from diffractive interac- 
tions. Differing from all other interaction mod- 
els GHEISHA creates the pions with momenta 
shifted to higher values for both displayed col- 
lision energies, i.e. with an excess in the range 
xiah > 0.4 and a deficit at xi a b < 0.3. In EAS sim- 
ulations this finding implies too high an elasticity 
which results in a stretching of the low-energy 
branches in the hadronic backbone. 

In p-air interactions we include the HSA [2U| as 
used in AIRES 2-6-0 [27| ■ It is parametrized only 
for air as target. Fig. shows the ^-multiplic- 
ities as function of the collision energy. Assum- 
ing FLUKA and UrQMD as standard we observe 
an overshoot of QGSJET as in Sect. O For 
GHEISHA the multiplicities seem now to be in 
agreement with the reference models. 
The HSA as implemented in AIRES has been 
" configured to approximately emulate the multi- 
plicities and energy distributions of other models^ 




Figure 6. Multiplicities of 7r + in p-air collisions 
as function of projectile energy. 
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Figure 7. Multiplicities of K + in p-air collisions 
as function of projectile energy. 



|27| . An interpolation between GHEISHA and 
QGSJET has been selected [22] to get a smooth 
transition to the high-energy model QGSJET. In 
my opinion HSA should not be tuned to imitate 
other models, rather the parameters should be 
adjusted in a manner to approach experimental 
values as closely as possible. 
The K + -multiplicities displayed in Fig. show a 
nice agreement of most models for energies be- 
low «50 GeV. In the HSA the K + -number has 
a significant deficit. By a comparison with the 
corresponding K~-multiplicities the discrepancy 
can be traced back to kaons resulting from the 
associated K + +A production which is missing in 
the HSA. GHEISHA exhibits a moderate deficit 
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Figure 9. Average transverse momentum of ti + 
in p-air collisions as function of projectile en- 
ergy. Experimental values of pp-collisions are 
taken from |25| . 



for both charged K-mesons which increases with 
the collision energy. Generally for negative kaon 
multiplicities the models show much better agree- 
ment than for the positive ones. 
Special care has to be taken of the transverse mo- 
menta as they might influence the lateral distri- 
bution of the EAS particles arriving at ground. 
In the HSA the transverse momentum distribu- 
tions need to be inserted from outside. For the 
other models we have examined the p± -feature al- 
ready in p- 9 Be collisions (see Fig.|3J|. Fig.[S]com- 
pares the p^-distributions of 7r + -mesons for all 
low-energy and several high-energy codes. The 
reference models FLUKA and UrQMD coincide 
up to p±< 1.75 GeV/c. Their predictions are well 
separated from that of HSA with its flat slope in- 



dicating an insufficient approximation to exper- 
imental values. This flat slope is detectable in 
the whole design range of the HSA within AIRES 
gEJ. The GHEISHA distribution exhibits even a 
completely unexpected peak at «1 GeV/c. As 
discussed in Sect. 13.21 with increasing p± the 
high-energy models fall-off with a too steep slope. 
The average transverse momentum {p±_) as func- 
tion of collision energy is given in Fig. Nearly 
all models follow the trend of experimental values 
of pp-collisions. HSA shows a too strong increase 
of (p_i_) with collision energy which results from 
the flat transverse momentum distribution shown 
in Fig.|Hl The parametrization of QGSJET below 
the lower end of its designed energy range results 
in too low an average transverse momentum, es- 
pecially below i?i a b < 500 GeV. 
For pion-induced collisions with air the transverse 
momentum distributions of protons arc displayed 
in Fig.^B Again the HSA shows a deficit at small 
transverse momenta and an overshoot at large 
px-values. The kink in the p^-distribution of 
GHEISHA at s=s 1 GeV/c once more demonstrates 
the problem in the kinematics of this model. 

4. AIR SHOWER SIMULATIONS 

The low-energy hadronic interaction codes have 
been combined with the high-energy code 
QGSJET 01 and in part with SIBYLL 2.1 to sim- 
ulate p-induced showers of 10 19 eV at vertical in- 
cidence. In the simulations the thinning method 
is employed with the parameters given in 
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Figure 11. Longitudinal profiles of charged par- 
ticle numbers for different combinations of inter- 
action models. 



Table CI For UrQMD the thinning level has been 
modified to e=10~ 5 because of the long required 
CPU-times. For all combinations the transition 
energy between high and low-energy model has 
been taken as 80 GeV. 100 showers of each model 
combination have been averaged to minimize the 
influence of shower-to-shower fluctuations. 

4.1. Longitudinal profiles 

The resulting longitudinal profiles are reproduced 
in Fig. [ill All averaged profiles with QGSJET 
reach the shower maximum at 777 ± 10 g/cm 2 . 
The differences exceed slightly the standard de- 
viation of 6.5 g/cm 2 . UrQMD with enhanced 
thinning fluctuations reveals the smallest value 
at 766 g/cm 2 , while FLUKA penetrates deep- 
est (X max = 786 g/cm 2 ). For the combinations 
with SIBYLL the value Y max = 791 g/cm 2 of 
FLUKA agrees well with 796 g/cm 2 of GHEISHA 
(for the influence of high-energy hadronic inter- 
action models on Y max see Fig. 6 in Ref. |18|L 
One might argue the differences coming from dif- 
ferent production cross sections, but for nucleon 
and pion projectiles with pi a b > 1 GeV/c the cor- 
responding cross sections agree within ±5 %, only 
the kaon cross sections show somewhat larger dif- 
ferences. No correlation of the X max - values with 
the cross sections of the low-energy models could 
be observed. 

4.2. Lateral distributions 

The lateral particle number distributions fall off 
by 5 orders of magnitude in the distance range 
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Figure 12. Lateral distribution ratios of par- 
ticle numbers for 7 (top), e (middle), and 
fj, (bottom) relative to the model combination 
QGSJET 01 + GHEISHA 600. 



from 200 m to 5 km. Therefore the ratios relative 
to the combination QGSJET + GHEISHA 600 
are displayed in Fig. 1121 This combination of 
high and low-energy models shows the flattest 
distribution for 7, e*, and /x^-densities. The 
reason might be the wrong kinematics resulting 
in a prolongated shower development at the low- 
energy end of the hadronic backbone. Close to the 
shower axis this effect is masked by the numerous 
particles emerging from high-energy interactions, 
while at large core distances with a longer shower 
development path and therefore larger hadronic 
shower age the excess in elasticity of GHEISHA 
600 produces more secondaries which finally flat- 
tens the lateral distributions. Similarly the high 
elasticity of GHEISHA 2002 found in Sect. El 
for 7T + - 14 N collisions may explain the flatter dis- 
tributions (with QGSJET and SIBYLL) for the 
/X -density. 
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Figure 13. 
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the lateral Cherenkov densities given in Fig. [T31 
are derived. For the Auger tanks of 
10 m 2 area the signal S(1000) of a verti- 
cal 10 19 eV shower (p-induced) amounts to 47 
(QGS JET+FLUK A) , 46 (QGSJET+HSA), 45 
(QGSJET+UrQMD), 43 (QGSJET+GHEISHA 
2002), 39 (SIBYLL+FLUKA), and 36 (SIBYLL+ 
GHEISHA 2002) VEM. Radial distributions of 
Cherenkov densities are plotted in Fig. [21 as ra- 
tios which again display the large dependence on 
the high-energy model at 1 km core distance. The 
differences increase at larger distances still rele- 
vant to measurements of the Auger detector. 
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Figure 14. Lateral distribution ratios of 
Cherenkov densities relative to the model com- 
bination QGSJET 01 + GHEISHA 600. 

While FLUKA and UrQMD show similar dis- 
tributions for the ^-densities (and to a minor 
degree for the e ± -densities), the flat distribu- 
tions of HSA have to be attributed to the insuffi- 
cient pj^-parametrization. The combinations with 
SIBYLL exhibit a significantly lower /^-density 
even at core distances as small as 1 km. 
For the energy calibration of the Auger exper- 
iment [H] the signal S(1000) of the Cherenkov 
water tanks at 1000 m core distance will be used. 
Here we investigate the dependence of S(1000) 
on the low-energy model. With the assumptions 
for the Cherenkov water tanks that 

• em-particles deposit their full energy, 

• muons contribute with their full energy, but 
at maximum with 240 MeV, and 

• 240 MeV deposited energy corresponds to 
1 vertical equivalent muon (VEM) 



5. IMPLICATIONS 
ING REMARKS 



AND CONCLUD- 



The CPU-time requirements to treat single colli- 
sions rsp. complete showers is important for EAS- 
simulations. In Table Q CPU-times are collected. 
The HSA-code as the fastest has an unsatisfac- 
tory performance as demonstrated in the preced- 
ing sections and lacks predictive power. A revi- 
sion of its parameters is needed. 
As shown in previous sections the kinematics of 
GHEISHA 2002 is still not satisfying. 
The microscopic UrQMD model gives reliable re- 
sults at the expense of extreme long CPU-times. 
In this model no parametrized cross sections are 
used, rather the projectile hits on a 'disk' cho- 
sen large enough to cover the maximum colli- 
sion parameter. Several trials (with their time- 
consuming microscopic calculation) without in- 
teraction are usually needed before in an inelastic 

Table 1 

CPU times for DEC-alpha 1000XP (500 MHz). 
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QGSJET 01, 6 = 0°, B h > 300McV, B M > lOOMcV, 
B cm > 250koV, W cm < 10 4 , W h < 10 2 . 

H.J. Dreschcr has accelerated UrQMD 1.3c by a factor 15. 
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collision secondary particles are produced. There- 
fore this model is less suited for EAS simulations. 

The FLUKA model shows the best performance 
and it is recommended for EAS-simulations de- 
spite the longer CPU-times. Unfortunately its 
source code is not yet publicly available, only ob- 
ject codes for 5 types of machine/operating sys- 
tems are distributed at present. 
Finalizing one can state that measurable shower 
parameters are predominantly influenced by the 
high-energy interaction programs |18j and only to 
some minor extend by the low-energy models. 
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